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The correct description of nondynamic correlation by electronic structure methods not belonging to the multireference 
family is a challenging issue. The transition of to D,\ /, symmetry in H4 molecule is among the most simple 
archetypal examples to illustrate the consequences of missing nondynamic correlation effects. The resurge of interest 
in density matrix functional methods has brought several new methods including the family of Piris Natural Orbital 
Functionals (PNOF). In this work we compare PNOF 5 and PNOF6, which include nondynamic electron correlation 
effects to some extent, with other standard ab initio methods in the H4 D.m, / Dih potential energy surface. Thus far, the 
wrongful behavior of single-reference methods at the D2h — Oi/ ( transition of H4 has been attributed to wrong account 
of nondynamic correlation effects, whereas in geminal-based approaches it has been assigned to a wrong coupling of 
spins and the localized nature of the orbitals. We will show that actually interpair nondynamic correlation is the key 
to a cusp-free qualitatively correct description of H4 PES. By introducing interpair nondynamic correlation, PNOF6 
is shown to avoid cusps and provide the correct smooth PES features at distances close to the equilibrium, total and 
local spin properties along with the correct electron delocalization, as reflected by natural orbitals and multicenter 
delocalization indices. 


I. INTRODUCTION 

The correct description of nondynamic correlation effects 
is a challenging task for electronic structure methods. In 
wave function approaches, a multireference ansatz is needed 
to properly account for these effects. The computational 
scaling cost of such methods limits their use to systems of 
moderate size. Within density functional theory (DFT) the 
proper inclusion of nondynamic correlation effects is an 
open problem^ In practice, a broken-symmetry calculation is 
usually performed producing wrong spin densities^ 

An alternative to both wave function and DFT methods is 
natural orbital functional theory (NOFT)^— In recent years, 
several functionals have been proposed by reconstruction 
of the two-particle reduced density matrix ( 2 -RDM) in 
terms of the one-particle reduced density matrix ( 1 -RDM )4 
In particular, within the family of Piris Natural Orbital 
Functionals (PNOF)^ PNOF 5 — and PNOF6^ are among the 
best candidates to treat nondynamic correlated systems. They 
describe properly the dissociation limit of several molecules, 
recovering the correct integer number of electrons on each 
fragment upon dissociation— ill Both PNOF 5 and PNOF6 
belong to the family of orbital-pairing approaches, but the 
former only includes intrapair electron correlation while 
in the latter electrons on different pairs are also correlated. 
The inclusion of interpair electron correlation in PNOF6 
allows a better description of correlation effects and it also 
removes the symmetry-breaking artifacts that are present in 
independent-pairs approaches such as PNOF 5 when treating 
delocalized systems.— 

The purpose of this manuscript is to analyze the effect of 
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interpair electron correlation on the treatment of nondynamic 
correlation by investigating the performance of PNOF 5 
and PNOF6 and several standard ab initio computational 
methods. To this end we will examine the D ih /D 2 h potential 
energy surface of the planar H4 model (hereafter, simply 
PES). 

H4 has been extensively used to test single-reference post- 
Hartree-Fock methodsfilHA and geminal-based theories, 22 2 5 
Hartree-Fock, MP 2 and MP 3 show a spurious cusp on the 
PES of H4 as the system evolves from D2h to D^h symmetry. 
The cusp is the maximum energy value along the symmetry 
transition. Conversely, traditional coupled cluster (CC) 
methods predict a cusp but this cusp is a local miminum in the 
D2h — t, transition. Recently, Bulik et al. have shown that 
an improvement of the description of correlated systems can 
be also achieved by removing terms in traditional CC theory 1 
Variational CC approaches also improve this wrong behavior 
of the traditional CC implementations , 1 2 i 13 i 19 however, most 
of these approaches revert the local minimum to a local maxi¬ 
mum but most of them do not avoid the presence of a spurious 
cusp. Geminal-based theories predict a (maximum) cusp 
at the square geometry. Jeszenszki et al l 2 have attributed 
this failure to an insufficient account of spin couplings and 
the localized character of the orbitals. By including triplet 
components in the geminals, the orbitals become delocalized 
and the characteristic cusp vanishes, but the resultant PES 
is not completely smooth and wave function becomes spin 
contaminated. The authors also examined the local spin? 4-26 
of the system using different geminal-based approaches. 
Jeszenszki et a/— have found that singlet-coupled geminals 
fail to describe correctly local spins at the O j /, geometry. The 
inclusion of triplet components improve the results but the 
local spin values are not smooth along the PES. 

Thus far, the wrongful behavior of single-reference meth¬ 
ods at the D2h — Dih transition of H4 has been ascribed to 
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a wrong account of nondynamic correlation effects ^' 12 ' 13 ' 19 
whereas in geminal-based approaches, the spurious (maxi¬ 
mum) cusp has been attributed to a wrong coupling of spins 
and the localized nature of the orbitals^ We will show that 
actually interpair nondynamic correlation is the key to quali¬ 
tative cusp-free correct description of H 4 PES. By introducing 
interpair nondynamic correlation, PNOF 6 is shown to avoid 
cusps and provide the correct smooth PES features, total and 
local spin properties along with the correct electron delocal¬ 
ization, as reflected by natural orbitals and multicenter delo¬ 
calization indices. 


II. THEORY 


where J pq = (pq\pq) and JC pq = (pq\qp ) are the direct and 
exchange integrals, respectively and C pq = ( pp\qq) is the 
exchange and time-inversion integral . 29 Matrices A and II are 
auxiliary matrices proposed^ to reconstruct the 2-RDM in 
terms of the NO occupancies. The diagonal elements of these 
matrices are A pp = n p and II PP = n p . The off-diagonal 
elements of A and II determine the different implementation 
of the PNOFt (i =1 6 ) series. In particular, PNOF5 and 

PNOF 6 differ on the treatment of the interaction between 
electrons on different pairs. 

In PNOF5, when orbitals p and q belong to the same sub¬ 
space fig, the off-diagonal elements of A and II are A pg = 

n q n p and 


A. PNOF5/PNOF6 

In this section we will briefly review the formulation of 
PNOF52 and PNOF6.i2 Both PNOF5 and PNOF 6 belong to 
the family of orbital-pairing methods, which divide the spatial 
orbital space into subspaces (a set of orbitals) that contain two 
electron each. These methods couple each orbital g below the 
Fermi level ( F = N/2, where N is the number of electrons 
of the system) with N c orbitals above it, being the sub¬ 
space containing orbital g and its coupled counterparts. The 
original formulations of both functionals were introduced for 
N c = 1 but subsequently extended versions (N c > 1) were 
reported^ 7 ’— The sum rule for the occupation numbers (n) is 
fulfilled for each of the N/2 subspaces 


^2 n P = 1 C 1 ) 

p efi s 

where p denotes a spatial natural orbital (NO) and n p its 
occupation number. 

The PNOF5 and PNOF 6 energy expressions for a singlet 
state system can be written as 

fi =E B s +EEE £ «- ( 2 ) 

9=1 pen f qen g 

The first term of Eq. (O corresponds to the sum of energies of 
F independent pairs with energy E g , namely. 


E g = ]T n p (2U pp + J pp )+ J2 E i \M ( 3 ) 

pefi 9 p,qeflg,p^q 

where T~L PP is the matrix element of the kinetic energy plus 
nuclear-electron attraction terms and J pp = (pp\pp) is the 
Coulomb interaction between two electrons with opposite 
spins at the spatial orbital p. The term E// 1 contains the inter¬ 
action energy between electrons in different spatial orbitals p, 
and q. 


E pq — ( n q n p ~ Aq p ) (2 Jpq — K, pq ) + II q p C pq (4) 


n p9 


-y/n^, p = g or q = g 
y/njn ^,, p,q > F, 
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respectively, and they vanish when p and q belong to different 
subspaces. Consequently, the second term of Eq. [2]becomes 

F 

EE E £“‘(PNOF5) = n q n p (2 J pq - K, pq ). ( 6 ) 
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The expression above indicates that the interaction between 
electrons in different pairs is treated at the mean-field level. 
Therefore, PNOF5 lacks correlation between electrons in dif¬ 
ferent pairs. In contrast, the PNOF 6 A pq and II p9 matrices 
(when p and q belong to different subspaces these matrices 
do not vanish) include terms that account for interpair elec¬ 
tron correlation. The off-diagonal elements A pi] and I I p(J in 
PNOF 6 read as 


A qv 

n 9 P 
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q < F,p > F 

S-y 

iL qp 

q > F,p < F 

— 9 9 

e °n q n v 

e~ s (n q n p )^ 

q > F,p > F 


where h p is the hole (1 — n p ) in the spatial orbital p and S, 
7 p , ,S' 7 , and II 7 are defined as 
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Recently, PNOF5 has been proved equivalent to an 
antisymmetrized product of strongly orthogonal geminals 
(APSG)4LiE Conversely, PNOF 6 is not related to geminal the¬ 
ories but it keeps the orbital-pairing scheme, Eq. [T] In this 
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work we have used the N c = 1 version of the functionals. 
That is, each orbital subspace contains two spatial orbitals 
and then only N spatial orbitals are correlated. In this sense, 
both functionals take into account most of the nondynamic 
correlation effects, but while PNOF5 includes only intrapair 
correlation, PNOF 6 incorporates also the interpair correlation, 
through A and fl matrices defined in Eq. [7](see Eq0i 


B. Local Spin And Electron Delocalization 

Local spins can be obtained by decomposing the expecta¬ 
tion value of the total spin square operator (S 2 ) into atomic 
or fragment contributions as 

{S 2 )=Y j {S 2 )a+ ]T(S 2 >ab, (9) 

A A=fiB 

where ( S 2 )a is the local spin on fragment A and ( S 2 )ab 
accounts for the coupling between spins on fragments A and 
B. Recently some of us have presented a general formulation 
of the local spin that fulfills a set of physical constrains . 24 ' 25 
For singlet systems, the formulation reads as 

{S 2 ) a = j (2Tr( 1 DS A ) - Tr^DS^D)) (10) 

+i £ - l ]T r ij]kl s A s A 

ijkl ijkl 

and 

{S 2 ) AB = iwS&sfj - 

ij kl ij kl 

where 1 D, F, and S A are the spinless 1-RDM, the spinless 
cumulant of the 2-RDM, and the fragment orbital overlap 
matrix * 24 The correct description of local spins has been 
recently put forward as a stringent condition to test natural- 
orbital based cumulant matrix (or 2-RDM) approximations^ 
and has been used to characterize and quantify the diradical 
and triradical character of molecules . 33 ' 34 In this work, we 
will use the local spin analysis to study the effect of the 
interpair electron correlation in PNOF5 and PNOF 6 on the 
spin coupling of electrons located at different atoms. 

The calculation of electron delocalization among different 
fragments can be performed through the NO-weighted over¬ 
lap multiplications involving the different fragments. This is 
commonly known as Giambiagi’s multicenter index-^ and its 
expression reads 2 & 

Iabcd = ^^ J n i n j n k r iiS i jSji,S)? l Si i ( 12 ) 

ijkl 

The quantity has been successfully used to account for several 
multicenter delocalization phenomena including multicenter 
bonding^ conjugation effects^ and aromaticity . 39 ’ 40 


III. COMPUTATIONAL DETAILS 

In this work we have computed the D|/,/D 2 /, PES of 
H 4 employing the following methods: Hartree-Fock (HF), 
CC singles and doubles (CCSD), CCSD with perturbative 
estimation of triple excitations (CCSD(T)), complete active 
space self-consistent field CASSCF (with a 4 electrons in 4 
orbitals active space), PNOF5, PNOF 6 and full configuration 
interaction (FCI). This benchmark data set includes methods 
that mostly include dynamic correlation effects (CCSD and 
CCSD(T)) or nondynamic correlation effects (CASSCF) and 
will be used as benchmark references to measure the amount 
of dynamic and nondynamic correlation effects included in 
PNOF5 and PNOF 6 . 

All calculations based on wave function methods have been 
performed with the Gaussian0?4f set of programs except those 
at the FCI level that were performed with a modified version 
of the program of Knowles and Handy4s22 NOF calculations 
have been carried out using DoNOF program. The matrix 
elements of the kinetic energy, the nuclear-electron attraction 
energies, and the one- and two-electron integrals needed to 
perform the PNOF calculations have been obtained from 
GAMESS4i-£ The correlation-consistent aug-cc-pVDZ^ 
basis set has been employed for all the calculations. The local 
spin analysis has been performed using DMN ^ 7 to compute 
the 2-RDM and APOS I-313-^ to calculate the local spins 
using the topological fuzzy Voronoi cells to define the atomic 

regions^ 


IV. RESULTS 

The PES of H 4 is characterized using two parameters, R 
and 9 (see Fig. Q}. The former, controls the distance be¬ 
tween each H atom and the center of mass while the latter 
measures the angle formed by two neighbor H atoms and the 
center of mass (see Fig. Q}. At 8 = 90°, the system possesses 
D 4 h symmetry and two configurations with symmetries a 2 b 2 u 
and a 2 b 2 u become degenerate. By modifying 9 one can con¬ 
trol the degree of symmetry distortion with respect to the D^, 
(8 = 90°) structure, thus modulating the multireference char¬ 
acter (and hence the nondynamic correlation) of the system. 
In this sense, the H 4 PES represents a challenging system for 
most electronic structure methods as it combines nondynamic 
correlation and dynamic correlation effects. 

The relative energies with respect to the minimum energy 
at 9 = 70° for each method of the H 4 model keeping R 
constant for different distances and modifying 8 are shown in 
Fig. [2] The system is symmetric at 9 = 90° and it is described 
by two degenerate configurations, which correspond to the 
minimum HF solutions at 9 < 90° and 9 > 90°, respectively. 
The FCI curve has an energy maximum at 9 = 90° and the 
energy curve is smooth along the entire range of angles. The 
energy needed to change from 8 = 70° to the D u, geometry 
decreases gradually as the radius R increases until the PES 
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FIG. 1. 0,,,/D,,, H 4 model. 


becomes considerably flat. The CASSCF curves show the 
right qualitatively features, i.e., a maximum at 9 = 90° and 
a smooth transition from 9 = 70° to 0 = 110°. However, 
due to missing dynamic correlation energy that becomes 
important at the 9 90° and 6 <C 90° regions, CASSCF 

relative values are downshifted to lower energies. 

At 9 = 90° two configurations become degenerate and 
the HF solution presents symmetry-breaking artifacts that 
result in a maximum cusp in the energy profile ^- 7 Therefore, 
it is only natural that most post-HF single-reference methods 
based on the RHF reference also fail to qualitatively describe 
this PES. Although at small R values CCSD and CCSD(T) 
mimic the FCI PES, as the radius R increases first CCSD(T) 
(at R = 0.80A) and then CCSD (at R = 1.20A) break down 
and show a cusp of the PES at 9 = 90°, which —unlike 
the HF cusp— is a local minimum with respect to 9. Since 
CASSCF with a (4,4) active space shows a qualitative right 
result and dynamic-correlation-including methods produce 
an artifact at 9 = 90°, one attributes this feature to the lack of 
nondynamic correlation effects. Consequently, at short values 
of R and for the 9 values considered, the CC results are in 
perfect agreement with FCI. 

PNOF5 —a nondynamic-correlation-including method— 
shows a maximum cusp at 9 = 90°, like VCC,— OQVCCD 
and OQVCCD(T)r^ and the lately introduced CCDO and 
CCSD0, which are single-reference CC variants that exclude 
certain excitationsr This result suggests that PNOF5 is miss¬ 
ing some nondynamic correlation and it is only this fraction 
of nondynamic correlation that is responsible for the spurious 
cusp. 

On the other hand, PNOF 6 which —at variance with 
PNOF5— includes interpair correlation, shows a smooth PES 
for R < 1.5 A, suggesting that only interpair nondynamic cor¬ 
relation is actually needed to obtain a cusp-free, qualitatively 
correct description of the H 4 PES at values close to the min¬ 
imum energy geometry. When R = 1.70A and 1 .90/1, the 
PNOF 6 solution is not perfectly smooth. This behavior is due 
to the crossing of two solutions of the PNOF 6 equations as 
can be seen in Fig. [3] In this graphic, the minimum PNOF 6 
solution is showed in solid lines. One can see the crossing of 
two solutions at 9 ~ 80°, 90°, and 100° for R = 1.70A and at 
9 ~ 70°, 90°, and 110° for R = 1.90A. At large R only one 
solution (labeled Sol. 2 in Fig. [3]) of the PNOF 6 equations is 
found, there is no longer a crossing and the PES smoothness 
is recovered, the shape of PNOF 6 and FCI relative energies 


TABLE I. Relative energies (kcal/mol) as the difference between the 
absolute energies at 6 = 90° and 9 — 70° for different values of 

r(A). 

Method R=0.80 R=1.00 R=1.20 R=1.40 R=1.60 R=1.80 


FCI 

68.75 

61.54 

48.92 

35.58 

23.79 

14.72 

HF 

99.15 

99.43 

93.38 

85.20 

76.59 

68.33 

CASSCF 

66.61 

58.19 

45.44 

32.68 

21.66 

13.27 

PNOF5 

87.63 

78.67 

61.88 

43.67 

27.94 

16.47 

PNOF 6 

74.19 

68.98 

57.74 

44.48 

30.68 

18.17 


being almost indistinguishable. 

Table Ugathers the relative energies at 9 — 70° with respect 
to the energy at 9 = 90°. For R = 0.8 A, R = 1.00 A, and 
R = 1.20 A PNOF 6 improves PNOF5 (as to compared to 
FCI) by 13.44, 9.69 and 4.14 kcal/mol, respectively. At larger 
values of R, PNOF5 improves over PNOF 6 but the difference 
between them does not exceed 3 kcal/mol. CASSCF results 
are closer to FCI than PNOF 6 for all the distances. The 
difference attains its maximum at R = 1.20 A, in which 
CASSCF is 12.29 kcal/mol closer to FCI than PNOF 6 . These 
deviations put forward the current limits of PNOF 6 to fully 
account for correlation effects. 

In table UJ we collect FCI, PNOF5, and PNOF 6 absolute 
energies for R = 0.80A, 1.20A, and 1.70A. PNOF5 energies 
are in all cases closer to FCI than PNOF 6 . This is due to the 
repulsive electron-electron interpair correlation energy term 
that is included in the PNOF 6 functional. PNOF 6 improves 
qualitatively the shape of the PES, provides good relative en¬ 
ergies at the price of higher absolute energies. 

APSG, which is the antisymmetric wavefunction behind 
PNOF5,— has been shown to also exhibit this spurious max¬ 
imum cusp at 9 = 90°.— The failure of APSG has been at¬ 
tributed to the localized nature of its orbitals and the wrong 
account of spin coupling. Szabados and coworkers 22 , have 
demonstrated that APSG using delocalized orbitals, which 
correspond to a solution of the ASPG equations, eliminates 
the cusp. In Fig. Q]we plot the orbitals that arise from PNOF 6 
and PNOF5 at R = 1.0 A and 0 = 90°. PNOF5 NO are 
localized on H — H bonds and each bonding orbital is cou¬ 
pled with its antibonding counterpart. At this value of 9, the 
same picture with the orbitals horizontally localized is equiv¬ 
alent. On the other hand, the PNOF 6 NO present the expected 
delocalized character and mimic the canonic orbitals obtained 
in a HF calculation. Importantly, both solutions showed in 
Fig. [3]for R = 1.70A and 1.90A present delocalized orbitals. 
Unlike PNOF5, PNOF 6 equations do not lead to a stationary 
solution that corresponds to a set of localized orbitals. 

The inclusion of interpair correlation also affects the 
occupation numbers of the corresponding NO (see table Hill). 
For small values of R at the CASSCF level, the a g orbital 
remains almost doubly occupied along the PES. The hz,, is 
doubly occupied for 9 <C 90° and there is a smooth transition 
from these structures to the 9 90° ones in which the doubly 

occupied orbital is the b^ u . At 9 = 90° both orbitals become 
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FIG. 2. Relative energies in kcal/mol with respect to the lowest energy found for each method 9 


70°, along the D^h/D^h PES of H 4 . 


TABLE II. FCI, PNOF5, and PNOF6 H 4 absolute energies in a.u. for different values of 6 and R 


R. = 0.80 A 


R = 1.20 A 


R = 1.70 A 


9 ° 

70 

72 

74 

76 

78 

80 

82 

84 

86 

88 

90 


FCI PNOF5 PNOF6 
-2.20639 -2.16625 -2.14177 
-2.19498 -2.15467 -2.12979 
-2.18305 -2.14247 -2.11728 
-2.17071 -2.12970 -2.10437 
-2.15805 -2.11639 -2.09113 
-2.14523 -2.10259 -2.07771 
-2.13244 -2.08832 -2.06428 
-2.12005 -2.07358 -2.05117 
-2.10881 -2.05839 -2.03902 
-2.10023 -2.04273 -2.02916 
-2.09683 -2.02660 -2.02354 


FCI PNOF5 PNOF6 


-2.14307 

-2.13184 

-2.12105 

-2.11077 

-2.10104 

-2.09195 

-2.08364 

-2.07635 

-2.07046 

-2.06652 

-2.06512 


-2.11972 

-2.10876 

-2.09810 

-2.08772 

-2.07760 

-2.06773 

-2.05809 

-2.04864 

-2.03936 

-2.03019 

- 2.02111 


-2.06829 

-2.05554 

-2.04323 

-2.03143 

- 2.02022 

-2.00970 

-2.00005 

-1.99151 

-1.98445 

-1.97926 

-1.97628 


FCI PNOF5 PNOF6 
-2.04310 -2.03505 -1.93271 
-2.03759 -2.02993 -1.92448 
-2.03271 -2.02533 -1.91722 
-2.02840 -2.02118 -1.91097 
-2.02462 -2.01744 -1.90572 
-2.02135 -2.01404 -1.90165 
-2.01855 -2.01094 -1.89992 
-2.01626 -2.00809 -1.89823 
-2.01451 -2.00543 -1.89659 
-2.01339 -2.00294 -1.89499 
-2.01300 -2.00055 -1.89423 


degenerate in terms of occupancies. The PNOF5 bonding 
orbitals are almost doubly occupied along the PES while the 
antibonding ones remain almost unoccupied. No degeneracy 
is observed in this case. By including interpair electron 
correlation, PNOF 6 NO and occupancies qualitatively mimic 


the CASSCF ones. It is worth noting that at 9 = 90° the 63 ,, 
and 6211 do not have exactly the same occupancy for most of 
the values of R shown in Fig. [HI] This might indicate that 
the interpair description is not fully recovered by PNOF 6 . 
The second solution shown in Fig. [3] as PNOF 6 (sol. 2), that 
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FIG. 3. Absolute PNOF 6 energies in a.u. for R = 1.70A (top) and R = 1.90A(bottom). PNOF 6 (Sol. 1 ) and PNOF 6 (Sol. 2 ) stand for the 
two solutions that show a crossing and PNOF 6 (min.) stand for the minim energy solution of each value of 9. 


TABLE III. CASSCF(4,4), PNOF5, and PNOF 6 NO occupation numbers at 6 = 90° for different values of R. 


R(A) 

Til 

n 2 

n 3 

77,4 



CASSCF 


0.80 

1.939 

1.000 

1.000 

0.061 

1.00 

1.882 

1.000 

1.000 

0.118 

1.20 

1.795 

1.000 

1.000 

0.205 

1.50 

1.604 

1.000 

1.000 

0.396 

1.70 

1.458 

1.000 

1.000 

0.542 

1.90 

1.327 

1.000 

1.000 

0.673 

20.00 

1.000 

1.000 

1.000 

1.000 



PNOF5 


0.80 

1.923 

1.921 

0.079 

0.077 

1.00 

1.835 

1.835 

0.165 

0.165 

1.20 

1.704 

1.704 0.296 

0.296 

1.50 

1.472 

1.471 

0.529 

0.528 

1.70 

1.335 

1.335 

0.665 

0.666 

1.90 

1.229 

1.229 

0.771 

0.771 

20.00 

1.000 

1.000 

1.000 

1.000 



PNOF 6 


0.80 

1.971 

1.185 

0.815 

0.029 

1.00 

1.942 

1.197 

0.803 

0.058 

1.20 

1.894 

1.191 

0.809 

0.106 

1.50 

1.771 

1.150 0.850 

0.230 

1.70 

1.645 

1.110 0.891 

0.355 

1.90 

1.495 

1.068 

0.932 

0.505 

20.00 

1.000 

1.000 

1.000 

1.000 


becomes the minimum energy solution for certain values of 
9 when R = 1.70/1 and 1.90/1 and is the minimum solution 
found for larger values of R, presents perfect degeneracy in 
terms of occupation numbers of the b 3u and b 2u orbitals for 
all values of 8 and R. 

The wrong coupling between spins located in diferent 
centers of the molecule is one of the causes for the failure 


of singlet-couplet geminal approaches to describe the H 4 
system. Jeszenszki et al. have used the local spin analysis to 
show that the inclusion of triplet components in geminals im¬ 
proves the APSG results but spin contamination appears when 
the triplet component in the geminal becomes important. The 
local spin value of one H atom of the H 4 system is shown in 
Rg. a As the system approaches the D^h symmetry, there 
is an increase of the diradical character of the system and 
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FIG. 4. PNOF6 (left), and PNOF5 (right) natural orbitals of H 4 for 
R= 1.0 A and 9 = 90° 



FIG. 5. Local spin values of one of the H atom of H 4 at R = 0.80 A 
with respect to angle 9. 


the local spin on atom H grows. PNOF5 cannot reproduce 
this trend and the local spin remains almost constant along 
the PES, while PNOF 6 local spin values in H 4 are in good 
agreement with the FCI results. 


Finally, let us examine the multicenter delocalization in the 
D 2 h to D.ih transition . Fig. [ 6 ] shows that PNOF 6 values 
closely follow the FCI ones and give a maximum electron de- 
localization in the D^h structures, whereas PNOF5 shows a 
rather constant profile, clearly indicating its inability to delo¬ 
calize the electron density along the H 4 skeleton. 



FIG. 6. Multicenter Giambiagi indices, Eq.QT] along the D 2 h — D^h 
transition for R = 0.80 °A performed with PNOF5, PNOF6 and FCI 
natural orbitals. 

V. CONCLUSIONS 

The PES of the planar H 4 model has been com¬ 

puted at several levels of theory. Single-reference methods 
show a spurious cusp at the D^h structure that thus far was 
attributed to nondynamic correlation. PNOF5 (which affords 
a correct description of molecular dissociation and other in¬ 
trapair nondynamic correlation effects) also shows a spurious 
cusp at D^h, whereas PNOF 6 provides a qualitatively correct 
description of this phenomenon. 

Since PNOF5 and PNOF 6 mainly differ from each other by 
the inclusion of interpair correlation, the factors responsible 
for the spurious description of the H 4 PES can be 

narrowed down to missing interpair nondynamic correlation 
effects. Indeed, the inclusion of interpair correlation in the 
pairing-orbital NOFT ansatz is key to recover the delocalized 
orbitals picture, remove the spurious cusp in the PES and 
properly account for the coupling between the spins located 
at different centers. On the other hand, inclusion of more 
terms to fully account for electron correlation seems to be 
needed to recover the smoothness of the curves at R = 1 .70/1 
and 1.90-4, to obtain quantitative results, and to recover the 
important correlation effects that separate PNOF 6 results 
from FCI. We hope that this study will shed light on the 
effect of interpair electron correlation and pave the way to 
the development of new electronic structure methods within 
NOFT or methods based on geminal expansion of the wave 
function. Research in this direction is underway in our 
laboratory. 
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